arXiv:1507.07599v2 [physics.soc-ph] 22 Nov 2015 


Effects of local and global network connectivity on synergistic epidemics 


David Broder-Rodgers 

Selwyn College and Cavendish Laboratory, University of Cambridge, Cambridge, raQ 

Francisco J. Perez-Reche 

Institute for Complex Systems and Mathematical Biology, SUPA, University of Aberdeen, Aberdeen, uM 

Sergei N. Taraskin 

St. Catharine’s College and Department of Chemistry, University of Cambridge, Cambridge, 

Epidemics in networks can be affected by cooperation in transmission of infection and also con¬ 
nectivity between nodes. An interplay between these two properties and their influence on epidemic 
spread are addressed in the paper. A particular type of cooperative effects (called synergy effects) is 
considered, where the transmission rate between a pair of nodes depends on the number of infected 
neighbours. The connectivity effects are studied by constructing networks of different topology, 
starting with lattices with only local connectivity and then with networks which have both local 
and global connectivity obtained by random bond-rewiring to nodes within certain distance. The 
susceptible-infected-removed epidemics were found to exhibit several interesting effects: (i) for epi¬ 
demics with strong constructive synergy spreading in networks with high local connectivity, the 
bond rewiring has a negative role on epidemic spread, i.e. it reduces invasion probability; (ii) in 
contrast, for epidemics with destructive or weak constructive synergy spreading on networks of 
arbitrary local connectivity, rewiring helps epidemics to spread; (iii) and, finally, rewiring always 
enhances the spread of epidemics, independent of synergy, if the local connectivity is low. 
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I. INTRODUCTION 

Dynamical processes on networks is a subject of broad interdisciplinary interest and intensive study 00 . In 
particular, network models provide a unique framework to describe a wide range of spreading processes including 
spread of infectious diseases, social phenomena or biological species Such models assume a graph representation 

of systems with nodes (vertices) that can be in several states specific to the spreading process (e.g. a host infected 
by a pathogen or a patch occupied by some species). The state of nodes can change due to interactions with other 
nodes connected by links (edges or bonds) to a recipient. For instance, such interactions may represent transmission 
of infection, opinions, behaviours or ecological migrations. 

The chances for a spreading phenomenon to affect a large number of nodes in a network, i.e. invade a network, 
depend on both the dynamics of interaction between nodes and the topology of the network 0&0. The simplest 
type of interaction is a pair-wise interaction when only two nodes are involved in transmission, e.g. an infected 
donor and a susceptible recipient. The process of transmission can be characterised by two parameters, the rate of 
transmission (or probability of transmission in discrete-time description) and time of interaction (time of existence 
of contact/link between the donor and recipient) which can be a random variable. If the transmission rate is a fixed 
constant parameter then the transmission of infection is a homogeneous Poisson process, i.e. simple transmission. 
However, in real situations the transmission can be a much more complex process. In particular, both the life-time of 
nodes in different states and the transmission rates can depend not only on characteristics of the donor-recipient pair 
but also on characteristics of other nodes. In other words, cooperative effects due to multiple-node interactions, which 
are called synergistic effects below, can affect the values of parameters characterising transmission. For example, the 
transmission rate can depend on the number of infected neighbours of a recipient. This number can change with time 
throughout the course of an epidemic and thus the transmission rate can change with time. These abrupt step-wise 
changes in transmission rate due to cooperative effects have significant and non-trivial effects on the spreading process. 

The role of the network topology in the ability of epidemics to invade the network can be of crucial importance. 
For example, assuming simple transmission, an invasion can be much easier in globally connected networks such as 
complete or random graphs than in lattices where the nodes are locally connected to their nearest neighbours in space. 
However, it is not clear a priori what the effects of cooperative or interfering synergistic phenomena in transmission 
would be on invasion in networks of different topology. For example, if the local connectivity in lattices is reduced 
by the rewiring of some bonds which enhances global connectivity, would this necessarily result in increased invasion 
ability for epidemics with synergistic effects? This and other related questions are addressed in the paper. Before 
going to the description of our model, we give a brief overview of existing models accounting for simple and complex 
transmission in spreading phenomena in networks of different topology. 

Models assuming simple transmission dynamics have provided good insight into some aspects of the interplay 
between the features of transmission and network topology. This is indeed the case for network models for epidemic 
spread, which often assume that the transmission of infection between a pair of donor-recipient nodes is independent 
of the rest of nodes connected to the pair EMU. A similar assumption was made in some models for spread of 
social phenomena 0 ] . These models predict that invasions are facilitated in networks with small local clustering 0 - 
0 and long links that can act as bridges for transmission of infection 0, 00. The size of invasions with 
simple transmission is therefore minimised in regular lattices (with relatively large clustering and no shortcuts) and 
maximised in random graphs (small clustering and many shortcuts). Small-world (SW) network topologies bridge 
the gap between lattices and random graphs by means of either random bond rewiring fig . 0, 10 ] or adding random 
shortcuts 0 , 0 ]. The size and chance of simple invasions increase with the probability of either rewiring (i.e. by 
reducing local and increasing global connectivity) or adding shortcuts (increasing global connectivity) 0 - tl8l I0 - [24| . 

Many social and biological systems involve complex transmission dynamics which are often characterised by syn¬ 
ergistic effects for donor-recipient pairs of nodes. These effects are not captured by simple epidemiological models 
but they can significantly change the dynamics of spreading processes. The knowledge and understanding of complex 
transmission dynamics on spreading and its interplay with the network topology is rather limited and is a topic of 
active research 0 0-0] . 

Synergistic effects can be either constructive or interfering (destructive). Constructive synergistic effects from the 
neighbourhood of a recipient node were explicitly observed in experiments on the spread of behaviour [0 ] and fungal 
invasion [0 ]. Social reinforcement was proposed as a key synergistic effect making invasions of social phenomena more 
likely and larger in clustered networks than in random graphs (i.e. opposite to the predictions obtained assuming simple 
transmission). This conclusion was supported by models involving social reinforcement from multiple neighbours [26i . 
0, HH, HH- The authors of Ref. j]3lj went a step further suggesting that these types of invasions are, in fact, optimal 
on SW networks rather than in fully clustered lattices. Interfering synergistic effects associated with, e.g. behavioural 
responses to epidemic spread |25l I34l - l36l ] or competition for resources [29j, can also play an important role in spreading 
dynamics. 

Constructive and interfering synergistic effects are often described separately. A recently developed model [29] 
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provides a flexible framework to study any degree of constructive or interfering synergy in any ty pe of network. With 
this model, it was shown that synergy affects the size, duration and foraging strategy of spreaders [2f|l32| and can even 
result in explosive invasions [l?} of epidemics with and without node removal and for the Maki-Thompson model [38j 
describing social phenomena. For regular lattice models, it was found that synergistic effects on invasion are enhanced 
by increasing local connectivity [32j ■ 

In this paper, we study the combined effect of local and long-range connectivity on synergistic spread. To this 
end, we use SW network models with rewiring which account for ubiquitous geographical constraints present in many 
social and biological systems [48]. We study a synergistic SIR process on such networks. The SIR model was originally 
formulated to investigate the spread of infection in populations where infected hosts either die or become permanently 
removed. In this model, the nodes can be in three states: susceptible (S), infected (I) or removed (R). For spread 
of social trends (e.g. opinion or rumour), similar states can be used to distinguish between ignorant individuals 
(analogous to S), individuals that spread the trend (analogous to I) and individuals that stopped spreading (R). We 
demonstrate that synergistic spread is str ong ly affected by the network topology. It is found that, in agreement with 
studies on social reinforcement |26 h 281 |3li I33| , systems with significant rewiring tend to be more resilient to invasion of 
epidemics with sufficiently constructive synergy. In contrast, interfering synergistic spread tends to be more invasive 
in rewired networks. We show, however, that these typical trends are very much affected by local connectivity. In 
particular, and in contrast to results in I26l-(28l. l3ll [33j . we show that rewiring systematically leads to larger invasions 
for weak local connectivity. We illustrate these and other effects with numerical simulations and analytical results for 
a simple model [H| based on an approximate mapping SIR synergistic spread to uncorrelated dynamical percolation 
(such mapping is exact in the absence of synergy [1C| |39|, |40|, l49jb 

The structure of the paper is the following. The model is introduced in Sec. Dj] and results of its numerical analysis 
are given in Sec. Mil The analytical results are presented and compared with the results of numerical simulations in 
Sec. IIVI The conclusions are given in Sec. |V] Some technical details are discussed in Apps. □ Dll □ and El 


II. MODEL 

Let us consider a network consisting of N nodes arranged on a regular two-dimensional lattice with each node 
connected to the same number of nearest neighbours, q > 2. In particular, we studied honeycomb (q = 3), square 
(q = 4) and triangular (q = 6) lattices in which the bonds connecting nearest neighbours can be rewired with 
probability 4> to a randomly chosen node under the constraint of no self- or double bonds. The probability of rewiring 
was assumed to be independent of the states of the nodes (cf. Refs. [HI, Ho})- Two types of models for bond rewiring 
were considered: (i) Spatial small-world (SSW) networks with rewiring to a random node within a finite distance, 
R € [I?inin,f?max], where -R m in > 0 and i? max are parameters of the model which are assumed to be independent of 
the lattice size (see Fig.QJa)); (ii) SW networks with rewiring to any random node within the system (see Fig. |ljb)) 
which is the limiting case of a spatial-SW if R m in = 1 and R ma .x(L ) ~ L —> oo. Here, L is the linear size of the system. 
In SSW networks, the rewiring is local which contrasts with SW networks where it is global. The bond rewiring was 
performed in the following way. Consider for concreteness a square lattice (see Fig.[T|). It can be constructed by an N 
times repeated translation of a node with two bonds attached to it (horizontal and vertical), over a distance a = 1 in 
both the horizontal and vertical direction producing N = L x L nodes on a square grid. These two bonds attached to 
each node are then rewired with probability cf> (per bond) to a random node within the range, [i? m i n , -Rmax], subject 
to no double bonds. Such a bond-rewiring algorithm similar to that used in Ref. flfil does not preserve the degree 
distribution, in contrast to degr ee-p reserving algorithms used in analysis of homogeneous SW networks [El . [ 5 IT , and 
it results in degree distribution J20|, which differs from that of a random graph (see App. [U]for more detail). 

The dynamics of the SIR process is defined by the life-time r = 1 of any node in the infected state and rate Aof 
stochastic transmission of infection from node i, infected at time tj, to an attached node j in the S-state. The SIR 
process with simple non-synergistic Poisson transmission is described by a constant transmission rate, \j = a. This 
model has been extensively studied on networks and it is well established that an SIR epidemic spreads (invades the 
network) if the transmission rate is greater than a critical value, a > a c , which marks the epidemic or, equivalently, 
invasion threshold iiEm. In the presence of synergy, the transmission rate of infection from an infected node, 
i, to a susceptible neighbour, j. is a piece-wise constant function. Step-wise changes occur after infection or removal 
events involving nodes i, j or any other in their neighbourhood. Such changes can be conveniently incorporated in 
numerical simulations using the event-driven continuous-time algorithm described in App. □ 

Our aim was to investigate the effects of synergistic transmission and rewiring (both varying probability and range 
of rewiring) on the invasion threshold. We focused on a specific type of synergy associated with the number of infected 
neighbours of a susceptible recipient node (referred to as r-synergy in [291] b In this case, the individual transmission of 
infection from node i to node j occurs with the rate A ij ( t ) which depends on the number rij ( t ) of infected neighbours 
of recipient node j excluding the attacker i [29l . [32]. The individual transmission between node i and j starts at time 
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FIG. 1: (Color online) (a) SSW network obtained from a square lattice with lattice spacing equal to unity and side 
L = 21 by means of rewiring with probability (f> = 0.03 within R £ [1, 5]. The rewired bonds are shown as dashed red 
lines between two nodes shown by red (gray) circles. The bonds going over the boundary reappear on the other side 
due to periodic boundary conditions, (b) SW network obtained from the square lattice by means of infinite-range 

bond rewiring with probability <j> = 0.03. 


ti when the node i was infected and it stops at t[, the time of infection of node j. not necessary by node i, or the time 
of recovery of node i, i.e. t[ = ti + r. The number rij(t ) of infected neighbours of j can vary with time in a step-wise 
manner for t £ [tj, t(). The time locations of the steps correspond to the stochastic infection and deterministic recovery 
events for the neighbours of node j excluding i and they depend on the history of the system at t < ij. 

In particular, we analysed invasions for the following three functional forms of A,,- (t) given as implicit functions of 
nj{t): 

(i) Exponential rate, 

A ij (t)=ae^, (1) 

(ii) the linear approximation to the exponential rate in Eq. m, 

Aij(i) = (a + af3rij(t))0( 1 + /3rij(t)) , (2) 

(iii) and a linear rate 2!)|. 32], 

A ij(t) = (a + f3'rij(t))9(a + . (3) 

Here, the Heaviside function takes the value 9(x) = 1 for x > 0 and is zero, otherwise. The expressions for these rates 
are valid for t £ For times outside this interval, A,;, = 0. 

The rate a in Eqs. ([T)l-© refers to the inherent (synergy-free) transmission rate. The coefficient /3 in Eqs. Jl]) 
and J2]) accounts for constructive (/? > 0) or destructive (/? < 0) synergy and is assumed to be independent of a. If 
[3 > 0 (/3 < 0) then the rate A ij(t) can exceed (be smaller than) the inherent rate a in the presence of a finite number 
of infected neighbours of node j (see Fig. [2]). The choice of the exponential function ensures positive values of the 
transmission rate for all values of /3. In the linear approximation to the exponential rate, the synergy contribution 
is proportional to the product of the inherent rate and synergy coefficient, i.e. oc a/3. However, it is possible that 
the synergistic effects do not depend on the inherent rate which is described by the functional form given by Eq. © 
where the synergy rate, /?', is independent of a. 


III. NUMERICAL RESULTS 

In order to study the effects of the model parameters (/3 or /3', i? m i n , Umax, and </>) on the invasion threshold, a Cl we 
numerically analysed the SIR process with synergistic rates defined by Eqs. ©-© on both SSW and SW networks 
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FIG. 2: (Color online) Schematic illustration of (a) a non-synergistic attack from infected node 1 (red or gray) to 
susceptible node 0 (blue or light gray) when all other nearest neighbours (nodes 2, 3, 4) of node 0 are susceptible 
(blue or light gray) and (b) synergistic attacks from nodes 1, 2 and 3, all infected (red or gray), to susceptible node 
0. The synergy effect in the case of simultaneous attack from three nodes in (b) is taken into account by a change in 
transmission rate from Aio = a in the case of a single non-synergistic attack in (a) to Aio = A 20 = A 30 = ae 2/S for the 

exponential form given by Eq. ©■ 


with periodic boundary conditions. All our simulations correspond to a linear size L < 200 for the underlying lattices. 
The transmission dynamics were modelled as a continuous-time synchronous process by using kinetic Monte-Carlo 
simulations as described in App. [A] 

For any given set ofparameters, the value of a c was estimated by using finite-size scaling analysis for one-dimensional 
spanning epidemics [29, in the case of SSW networks with finite-range rewiring, R max <C L. For SW networks 

with infinite-range rewirin g, w e used a linear fit [ 21 j and/or inflexion point for the mass of the infinite cluster vs 
inherent transmission rate 111 Il7j (see App. [B] for details). 




P P 


FIG. 3: (Color online) Invasion phase diagrams (inherent critical rate, u c , vs synergy parameter, /3) for SIR 
epidemics on (a) a SSW built from a square lattice (q = 4) with rewiring range R € [1,13], 

(b) a SW network with infinite-range rewiring on a square lattice. For both finite- and infinite-range rewiring, the 
exponential form of the transmission rate given by Eq. © was used. Different line styles refer to the phase 
boundaries corresponding to different values of 4> as marked in the figure legend. The crossing points of the phase 
boundaries for the two limiting cases of complete (solid lines for <fi = 1 ) and no (dashed lines for </> = 0 ) rewiring 
occur at /3 = /?** where /3** ~ 2.9 ± 0.3 in (a) and /?** ~ 2.5 ± 0.5 in (b). 


Fig. [3] shows how the critical value of the inherent transmission rate depends on the synergy parameter (3 for the 
exponential rate given by Eq. © on a square lattice (q = 4). Results are shown for SSW (panel (a)) and SW (panel 
(b)) networks. Each line in the figures gives the invasion threshold a c as a function of (3 for given </> and rewiring 
range. SIR epidemics in systems with a and (3 above/below the invasion line are invasive/non-invasive. For fixed 
values of (/, the critical transmission rate A c depends on two parameters a c and f3 in such a way (see Eqs. ©-©) that 
if /3 increases, then the value of a c should decrease in order to keep the same value of A c . Therefore, as expected, a c 
decreases with increasing (3 for any fixed </><!, meaning that increasing synergistic cooperation systematically makes 
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FIG. 4: (Color online) The dependence of the critical inherent transmission rate on rewiring probability <f> for 
different values of /3 between /3 = 5 (lowest curve) and /3 = —2 (uppermost curve) with the other curves 
corresponding to values of /3 incremented by unity (the lines are shown to guide the eye only). Panels (a) and (b) 
show data for the same models as in Fig. [3] The horizontal lines correspond to a c ~ a* and /3 ~ /3* with 
(/3*, a*) ~ (1.88 ± 0.04,0.46 ± 0.02) in panel (a) and (/?*, a*) ~ (2.03 ± 0.03,0.40 ± 0.03) in panel (b). 


systems less resilient to epidemic invasion. The effect of rewiring probability on a c is more involved and the trend 
depends on the value of (3. For small rewiring probability (<)> <C 1), all the phase-separation lines intersect at a single 
model-dependent point (/3», a*). This means that for a certain value of /3 = /3*, the critical inherent transmission rate 
a c = a* does not depend on the rewiring probability </> (see horizontal lines in Fig.[]]and dependence of the mass of the 
infinite cluster on a in Fig. [TOl' bl in App. [Bj . For values of /3 < /3*, the critical threshold in a decays with increasing 
(f> (see Fig. 0] and Fig. flOl all. This is the expected behaviour for synergy-free epidemics in SW networks [Hj]. In 
contrast, for relatively strong synergy, /3 > /3*, the critical value of a increases with rewiring probability (see Fig. 0] 
and Fig. HOf clh In other words, the more bonds rewired in the system, the more resilient it becomes to strongly 
synergistic SIR epidemics. This finding is in qualitative agreement with the social reinforcement effects discussed 
by Centola [27j |. It is important to note, however, that this effect does not only require the synergistic effects to be 
constructive (i.e. /? > 0) but it also requires that they are strong enough so that /3 > > 0. 
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FIG. 5: (Color online) The dependence of the critical inherent transmission rate on the probabilities (a) </> cu t of just 
cutting local bonds and (b) </> a dd of just adding global bonds for different values of /3 on a square lattice (same line 
styles as in Fig. 01 ). The exponential form of the transmission rate given by Eq. © was used for obtaining the data 
presented in both panels and bonds were added in the finite range R £ [1,13] for the model in (b). 


In order to interpret the results presented above, it is helpful to analyse separately the effects of cutting local bonds 
and adding local or global connections which are the two basic operations involved in rewiring. For large pos itive 
/ 3 , removing short-range bonds alone increases the resilience of the system whilst adding bonds on its own 11, [2lJ 
decreases its resilience to SIR epidemics (see Fig. 0]). In the rewiring scenario, these two tendencies compete and the 
resulting effects on invasion depend on the synergy strength. The effect of enhanced resilience with increasing rewiring 
is observed only for relatively large values of the synergy parameter /3 when the constructive synergy helps the SIR 
process to evolve locally where support from the infected neighbours is strongest. The rewired bonds diminish the 
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local connectivity and bring the infection to such remote parts of the system where there are practically no infected 
nodes which could support further spread of infection. The loss in ability to spread locally is more significant than 
the gain due to jumps to remote places where the advantages of high constructive synergy cannot be used efficiently. 
This is the reason why addition of new rewired shortcuts can make the system more resilient. 

The arguments presented above relied on the interesting property of the intersection point (0*,a*) which does not 
significantly depend on rewiring with 0 ■Cl. However, the conclusions remain valid for any value of 0. In general, the 
intersection point of the phase boundaries for arbitrary 0 £ (0,1) and for 0 = 0 occurs at a point (0(0), d(0)) which 
depends on 0. The function d(0) monotonically decreases with 0 from d(0) = a* to d(l) = a** < a*. In contrast, 
0 (0) is a monotonically increasing function taking values between 0(0) = 0* and 0(1) = 0** > 0*. These trends can 
be seen in Figs. EJa) and[6])a) for networks with underlying square and triangular lattices, respectively. 

Considering a fully rewired lattice with 0 = 1 allows the effect of enhanced resilience to SIR epidemics with strong 
synergy to be predicted independently from the analysis of SW networks. Indeed, for 0 = 1, the network becomes 
similar to a random (Erdos-Renyi) graph where local lattice connections are rare, resulting in an absence of small 
loops. Therefore, in this limit, the simultaneous presence of more than one infected node in the neighbourhood of the 
recipient is very unlikely and thus the synergy effects should be negligible. In other words, the epidemic threshold 
does not depend significantly on 0. Consequently, the phase boundary is close to a horizontal line in the 0-a plane 
and, importantly, this line can intersect the phase boundary for 0 = 0 at the point (0**, a**) (see Fig.[3]and Fig.[HKa)). 
This means that for 0 > 0„* the fully-rewired system is more resilient to invasion than the original lattice without 
rewiring. Again, this result contrasts with the prediction from models with simple transmission suggest ing that, given 
a mean node degree, invasions are more likely in random graphs than in regular lattices ll6Hl8i [21H24} . For instance, 
the critical transmissibility, T c = 1 — e~ a ° (a complementary quantity marking epidemic threshold, see Sec. IIVI for 
more detail), on a square lattice ( q = 4) T c = 1/2 whereas it is T c = 1 /{q — 1) = 1/3 in a random graph with 
( k) = q = 4. Our model reproduces this traditional behaviour for 0 < 0**. Such an effect becomes more pronounced 
for regular lattices with higher coordination number. For example, in a triangular lattice (q = 6), all the epidemics 
with exponential transmission rate are invasive for relatively large values of 0 > 0 m ax — 4. However, in the case 
of a fully-rewired triangular lattice, the presence of an almost horizontal phase boundary at a c = a ~ a** (with 
a** ~ 0.21 ± 0.03 for both the SW network and SSW network with R £ [1,13]) clearly demonstrates that all the 
epidemics become non-invasive for a < a** including those which were invasive in non-rewired networks for 0 > 0 max 
(see Fig. Oa)). 




0 0 

FIG. 6: (Color online) Invasion phase diagrams for (a) triangular (q = 6) and (b) honeycomb (q = 3) lattices with 
finite-range rewiring with R £ [1,13]. An exponential rate given by Eq. (|TJ) was used. The same line styles as in 
panel (a) for 0 £ [0,1] were used in panel (b). The crossing point in (a) between the phase boundaries for networks 
without rewiring (0 = 0) and complete rewiring (0 = 1) corresponds to 0 = 0„* ~ 0.9 ± 0.1. 


In lattices with relatively small coordination number (i.e with weak local connectivity), the effect of the rewiring- 
enhanced resilience to synergistic epidemics is not observed. This is due to the fact that the effects of synergy in 
such lattices are not very strong and the phase boundaries in the lattices without rewiring are almost horizontal. 
This effect is illustrated in Fig. EOT) for the honeycomb lattice (q = 3), where the invasion threshold is in the range, 
1.02 < a c < 1.15. In the fully rewired lattice, the critical inherent rate takes values a c = a ~ 0.57 ± 0.07 for infinite- 
range rewiring and a c = a ~ 0.62 ± 0.06 for finite range rewiring with R £ [1,13]. In both cases, a c is practically 
independent of 0 and the invasion boundary is an almost horizontal line located at d, below the range of the invasion 
boundary corresponding to 0 = 0. Therefore, it is not surprising that rewiring decreases the resilience of the system 
irrespective of the value of 0 (see Fig. EDb) ). These results show that the effect of local social reinforcement in networks 
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with weak local connectivity is not strong enough to compete with the gain in the spread efficiency achieved by the 
shortcuts and rewiring in the networks. This effect is not captured by existing models with social reinforcement [26}- 

0,0, ip. 

The critical value of d, defined by the position of the horizontal phase boundary in the (/3, a) plane for fully rewired 
lattices, can be found in terms of the bond-percolation threshold |l0|, T c , as 

d = ln(l^T c )- 1 , (4) 

where T c = (k)/(k(k — 1)} [56|-[58[. The values of a obtained numerically for infinite-range fully rewired lattices with 
triangular (d ~ 0.21 ± 0.03), square (d ~ 0.35 ± 0.04), and honeycomb (d ~ 0.57 ± 0.07) geometries agree well with 
the values calculated from Eq. © of d ~ 0.201 ± 0.001, d ~ 0.336 ± 0.002 and d ~ 0.558 ± 0.005 for the same lattice 
types, respectively. 

The effect of rewiring-enhanced resilience to invasion does not change qualitatively for the variety of the models 
given sufficiently high coordination number of underlying lattice. In particular, we found it for different ranges of 
rewiring including relatively small ones. Similarly, the effect was observed for models with linear transmission rates 
given by Eqs. m and ©. Fig. [7] summarises our findings for the variety of models investigated. Here, we show the 
location of three sets of characteristic points (/3*,a*) found for different models with exponential (Eq. (JT]) ), linear 
approximation to the exponential (Eq. ©) and linear (Eq. ©) transmission rates defined on square lattice. These 
points, as expected (see Sec. lIVI) . belong to the corresponding phase separation lines, a c (/3, $ = 0), for models without 
rewiring (<^> = 0). 

It should be noticed that the position of the characteristic point (/3*,a*) changes in a systematic way, moving 
down along the phase boundary (for lattices without rewiring, see solid curve in Fig. © with increasing rewiring 
range, tending to the point corresponding to infinite rewiring range. This can be understood as follows. For given 
/3 and (f> <C 1, increasing the maximum rewiring radius, I? max , will increase the global connectivity whilst the local 
connectivity remains similar. This has the effect of making the system more susceptible to invasion and thus, the 
larger the rewiring range, the smaller the critical inherent rate. Therefore, the crossing point with the monotonically 
decaying curve corresponding to a c (/3,4> = 0) shifts downward along this curve. 



FIG. 7: (Color online). The set of characteristic parameters (/?*,a*) for several SSW models and the SW model 
(with small probability of rewiring, 4> < 4.5 x 10~ 2 <C 1) defined on a square lattice with (a) exponential (Eq. (©)), 
(b) linear approximation to the exponential (Eq. I©) and (c) linear (Eq. ©) transmission rates of infection for SIR 
epidemics. In each panel, the solid line shows the invasion threshold separating non-invasive (below the line) and 
invasive (above the line) regimes in networks without rewiring (</> = 0). In panel (a), the continuous line corresponds 
to part of the curve in Fig. [3] Different symbols refer to different models of rewiring within the range 

R £ [7?min, 7? max ] as marked in the figure legend. 


IV. ANALYTICAL RESULTS FOR A MODEL WITHOUT CORRELATIONS IN TRANSMISSION 


Synergistic epidemics with removal can be viewed as correlated dynamical bond-percolation [29|. In this mapping, 
the bond probability between two nodes i and j corresponds to the probability that node i infects j during its infectious 
period, r, 





Tij = 1 - exp 


(5) 
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where ti is the infection time of node i. The probability is usually called the transmissibility from node i 
to node j. Synergy makes the transmission rate \j{t ) dependent on the infection history of the pair i — j and its 
neighbourhood. Accordingly, transmissibility will be in general different for different pairs of hosts, i.e. transmissibility 
is heterogeneous over the set of pairs of hosts. This heterogeneity is annealed, i.e. it varies with time, as opposed 
to quenched heterogeneity which is well studied for epidemics |40l - l47j . and its effect on the spread of epidemics is 
not obvious. In addition, transmissibilities of sufficiently close donor-recipient pairs have common nodes in their 
neighbourhoods and are not independent from each other if transmission is synergistic (since the infection history of 
close donor-recipient pairs neighbourhoods can overlap). In Ref. [29j . we showed that correlations in transmissibility 
can have a significant effect on invasion for large synergy. In spite of that, we found that the main features of invasion 
phase diagrams on lattices can be qualitatively described by a model which ignores spatial correlations but accounts 
for crucial spatial heterogeneity in transmissibilities [29l 1321 ]. Here, we extend this approach to obtain approximate 
analytical results that explain the rewiring-enhanced resilience reported above. 

The critical transmissibility in rewired networks, T c (<()), (i.e. the bond-percolation threshold) coincides with the 
mean transmissibility, (T((/>)), in the system: 

TO)> = T cW • (6) 

The value of T c (0) depends on the topology of the network through the rewiring probability, but it does not depend 
on a or 0 (0')- The expression for the mean transmissibility, (T((/>)), involves averaging over degree distribution, 

TO)) = I>< r *>, (?) 

k 

and averaging over possible challenge histories of recipients with fixed number k of nearest neighbours, accounted for 
by (Tfc) in Eq. 0. 

Eq. 0 is valid for heterogeneous transmissibilities [U, |4l|, [ssj but it assumes the absence of correlations in trans¬ 
missibilities for different bonds, which is true for non-synergistic SIR processes with a fixed removal time [43|. As 
argued above, such correlations are inherent for spread of the synergistic SIR process and condition © does not hold 
in general HqI, (H, EE EH ■ However, assuming that Eq. © holds even for synergistic SIR processes leads to a quan¬ 
titatively correct invasion phase diagram for small values of 0 (/?') and a qualitatively correct picture for relatively 
large values of 0 ~ 1 (0' ~ a) [29|, |32[ • In order to analytically study the consequences of Eq. (0 , we linearise the 
dependence of ( T(<j >)) and T c (</>) on </>. This leads to an approximate condition for epidemic threshold which reads as 
(see Ann. ID for more detail), 


T c0 (q) -M= 1 - e_ “ (! - » q {a, 0)H(/3')) 9-1 

- 2cj)e~ a [(1 - s q -i(a, 0', O)B(0')) g ~ 2 + (1 - s, +1 (a, /?', 0)R(/?')) 9 ] 
+ ^e- a {l-s q {a,0' 1 Q)B{0')) q - 1 

+ (9 - 1 )< Pe - a B ( 0 ') ds ^^ O) (1 - s ,( a , 0 ’, 0 ) B ( 0')) q ~' 2 , 


( 8 ) 


where T c o(q) is the bond-percolation threshold for a regular lattice with coordination number q and the non-negative 
functions s q (a, 0’, </>), B{0’) and A q are introduced in Aon. ICl 
Eq. 0 can be solved for 0') resulting in, 


T c o{q) 

<t>= - 


1 — e _ “ (1 — s q (a, 0', O)B(0')) q ~ 1 
' F q (a,0') 1 


(9) 


if F q (a, 0') 
F q (a,0 0 


^ 0, where 
= Aq- e~ a 


2(1 - s g _ 1 (a, 0', 0 )B{0')y~ 2 + 2(1 - s q+1 (a, 0', 0 )B{0')Y 


+ (q-l)B(0')(l-s q (a,0',O)B(0')y- 2 


ds q (a,0\ 0) 
d(l> 


4(1 -s q {a,p',n)B{P')) q - 1 


( 10 ) 


For a fixed value of 0', Eq. © defines how the critical value of the inherent transmission rate varies with rewiring 
probability. 
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For synergy-free epidemics with /?' = 0, the values of s q = 0 and F q (a, 0' = 0) = A q . Accordingly, the critical 
inherent rate decreases with rewiring probability, 


Oic — cUo In 



Aq<t) \ 

1 -T c0 J 


ot c o 


1 - T c0 ’ 


( 11 ) 


where a c o = — ln(l — T c q) is the critical transmission rate in the lattice without rewiring and synergy. However, 
for increasing synergy, the decay of a c with increasing rewiring probability becomes less pronounced and eventually, 
depending on Sk(a, /?'), it can become an increasing function. This happens at a characteristic value of /3' = /?', when 
a c = a* does not depend on <j>. Within the linear approximation for <f> -C 1, this is possible when both the numerator 
and denominator in Eq. © are simultaneously equal to zero, i.e. 


T c o — 


1 — e 


(l-s 9 (a*,/3;,0)H(/3:))' 3 


= 0 


F q (a*,/3i) = 0 . 


(12a) 

(12b) 


On the one hand, Eq. (U2al) gives the phase boundary for invasion in the absence of rewiring. On the other hand, the 
condition imposed by Eq. (I12bl) ensures that the solution (/3*,a») does not depend on (f>, as was found numerically 
(see Fig. 0. 

For given /?', the solution of Eq. © for a c {4>) qualitatively agrees with numerical simulations. In the case of 
triangular (see Fig. ©a)) and square (see Fig. ©b)) lattices, the change from a decrease of ot c with increasing <j> for 
/3' < /?* to an increase for /?' > /3' is observed. At a characteristic value /3' = /?', the inherent rate does not depend 
on (f>. For honeycomb lattices (see Fig. ©c)), no such transition is seen and the critical inherent rate only decreases 
with rewiring probability (j>. 

The model is accurate for /? = 0 when there are no synergy effects and thus, no correlations in the transmission 
between different pairs of nodes (see top lines in Fig.©. Significant deviations between numerical data and model 
predictions are seen for larger values of ft > 1 and they are due to approximations ignoring correlations in transmission. 
In spite of that, the model still provides qualitatively correct tendencies in a c (</>) for different values of ft' in various 
lattices. 

Interestingly, in the case of just cutting or adding bonds, the condition given by Eq. (112bl) is not satisfied for any fj', 
meaning that a c varies with <f> as a monotonically decreasing or increasing function for just cutting or adding bonds, 
respectively (see Ann. iDl). 


V. CONCLUSIONS 

To conclude, we have investigated the effects of local and global connectivity on spread of synergistic epidemics. 
The underlying networks used in the analysis were two-dimensional lattices with different coordination number (hon¬ 
eycomb, square and triangular). The local and global connectivity in these networks were changed by means of local 
(finite-range) and global (infinite-range) random bond rewiring. The global bond rewiring produced two-dimensional 
small-world networks while the local rewiring created spatial small-world networks with geographical constraints 
on the finite length of rewired bonds. SIR epidemics with constructive and destructive synergy transmissions were 
analysed on such networks. Our main findings are the following: 

(i) Bond rewiring enhances resilience to synergistic epidemics if two conditions are satisfied. First, the synergy 
effects are sufficiently strong and, second, the local connectivity is high enough. More specifically, the effect 
of rewiring-enhanced resilience is found only on lattices with high coordination number (q > 4) and synergy 
strength (3 > /?* >0. This finding is in line with those in Refs. (2614281 (31. 33]. 

(ii) Independent of local connectivity, if constructive synergy is not strong enough, i.e. 0 < /3 < /3*, or synergy 
is destructive (/3 < 0), rewiring enhances the spread of (reduces the resilience to) epidemics. In other words, 
destructive and weakly constructive syner gy d o not change qualitatively behaviour of synergy-free epidemics 
in rewired (small-world) networks (l7j, [u| |21, |2j, 124J. In particular, the fact that the traditional framework 
is recovered for 0 < /3 < /3* challenges the statement of Refs. (2614281 l3li . [33 ] showing that rewiring-enhanced 
resilience of epidemic invasion does not occur for every constructive synergistic mechanism. 

(iii) Independent of the strength of the synergy (constructive or destructive), if the local connectivity is small enough, 
the rewiring always decreases the resilience of the network to SIR epidemics. In particular, we have demonstrated 
this effect for epidemics in rewired honeycomb lattices (q = 3). 
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(j) (f) 



<t> 


FIG. 8: (Color online) The dependence of the critical inherent transmission rate on the rewiring probability 4> for (a) 
triangular (q = 6), (b) square ( q = 4) and (c) honeycomb ( q = 3) lattices. The results were obtained for a linear 
transmission rate given by Eq. (© and rewiring within a finite range of RG [1,4]. The symbols refer to the 
simulation results whilst the lines represent the model results given by Eq. ©. The same symbols and line styles are 
used in panels (b) and (c). The data for only two values of /3' are shown in (a) because a c —> 0 for /3' > 0.65 (see 

Ref. [H). 


All these three effects are quite robust to changes in the functional form of synergy transmission rate. In our 
approach, synergy is modelled by the dependence of the transmission rate of infection between a donor-recipient pair 
on the number of infected neighbours of the recipient. Three types of functional dependence of the transmission rate 
on the number of infected neighbours of donor-recipient pairs were investigated. Similar effects of local and global 
connectivity on spread of synergistic epidemics were found for all of them. The strength of synergistic effects was 
controlled with a single parameter, /3, which allowed both constructive (/3 > 0) and destructive (/3 < 0) synergy 
effects to be studied. This might be considered as an advantage of our model relative to other approaches typically 
studying one type of synergy, either constructive or destructive [H, iH, |3jJ. The effects reported here correspond to 
small-world networks obtained with a rewiring strategy which brings heterogeneity in the node degree. However, such 
heterogeneity is expected to play a secondary role on synergistic effects (e.g. rewiring-enhanced resilience) compared to 
rewiring-induced changes in local and global connectivity. In particular, we expect similar interplay between synergy 
and local/global topology when using a rewiring strategy leading to small-world networks with homogeneous node 
degree (22M5l|. 


ACKNOWLEDGMENTS 

We would like to thank Will Jennings for participating in developing the code used in simulations. 


Appendix A: Algorithm 

In this Appendix, we describe the rules of the SIR process and algorithm used in the simulations. 
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The SIR process can be described as a trajectory in discrete state space with the state vector having N components 
which can have three discrete values corresponding to different states of the nodes (S, I, and R) resulting in 3 N states 
in total, {5i}. The process evolves by means of instantaneous jumps between the states Si. These jumps occur at 
times ti (elapsed from the start of the process at t = 0) when the system rests in state Si and the trajectory is an 
ordered in time sequence of states <Sj(ij). The jumps can occur only between the states described by state vectors 
which differ in one component only, i.e. only one node changes its individual state after the jump. Only one of two 
changes are possible in one event: infection, i.e. S—>-1 or removal, i.e. I—>R. 

The time intervals between jumps are defined by the dynamical rules of the process, namely by the rules for infection 
and removal. The removal rule states that a node becomes deterministically removed after time r (parameter of the 
model) elapsed since the moment of infection of this node. For example, if node j has been infected at time t t with 
the system being in state Si , it is removed at time tk =U + t when the system is in the state Sk ■ The states Si and 
Sk can be separated on the trajectory of the process by many other system states which correspond to infection and 
removal of other nodes. The individual transmission of infection from an infected node to a susceptible one connected 
to the infected node occurs stochastically at times given by a Poisson process. The rate of such processes remains 
constant during quiescent intervals of time between consecutive states but, in the presence of synergy, its value can 
change after transition events. For example, assume that the system jumps from the state <Sj_i to the next state on 
the trajectory Si at time ti and infection event occurs at time ti , i.e. a susceptible node k becomes infected. The 
process of infection of node k is a superposition of independent individual infection transmissions from all infected 
neighbours connected to k and it takes place with the rate A k = X), n w h ere m runs over all infected neighbours of 
k. The value of A k does not depend on time for t £ [ti- 1 , ti) and it is fully defined by the state of the system <Sj_i at 
t = ti- 1 , i.e. it does not depend on the previous history of the system at t < ti- 1 . In particular, the rate A k depends 
on the number of infected neighbours and on individual rates X m k- In general, the rates A m k can also depend on the 
number of infected neighbours of node k (only for non-synergistic epidemics the values of A m k are independent of the 
infected neighbours of k). After the infection event at time t = ti, infection rates between any infected-susceptible 
pairs of connected nodes may have changed and should be updated. Similarly, if a deterministic removal rather than 
infection event takes place at ti then all the individual infection rates should also be updated. 

Numerically, we aimed to sample without bias all possible trajectories of the SIR process. This can be achieved 
by means of kinetic Monte-Carlo [59j simulations exploiting the Gillespie algorithm (direct method 00) with 
modifications accounting for deterministic recovery events. Within this algorithm the SIR trajectory was sampled as 
follows. 

1. Start simulations at t = 0 by infecting a small number of nodes, Nq ~ 0(1) -C N, distributed randomly within 
the network. Create a list of infection events, i.e. the list of susceptible nodes linked to the infected nodes and 
cumulative infection rate (sum of all individual infection rates) for each node in the list. Create a list of nodes 
in the infected state with their recovery times. 

2. For a current time step t , calculate the cumulative infection rate, R = )TA . Ay where i runs over all infected 
nodes and j runs over susceptible neighbours of infected nodes connected by the links to them, i.e. the sum is 
evaluated over all possible individual infection transmissions in the network. 

3. If R > 0, calculate a uniformly distributed random number r\ £ (0,1]. 

4. Calculate a time step till the next possible infection event, At = — ln(ri)/R. 

5. Compare t + At with the time of the earliest deterministic recovery event, t r (> t). If t + At > t r or R = 0 
perform the recovery event at t T , update the list of the individual infection rates, set the current time t to t = t T 
and return to step 2. 

6. If t + At < t r , calculate a uniform random number r 2 £ (0,1]. 

7. Add the individual infection rates for nodes from the list of infection events cumulatively until it exceeds r 2 -R. 
Infect the node for this event, update the infection rates and recovery times and set t = t + At. 

8. Return to step 2. 

The simulation stops when there are no nodes in infected state. This algorithm is valid for both synergistic and 
synergy-free SIR processes. For a synergy-free process, the infection rates Ay = a are identical for all the individual 
infection processes. In case of synergy, the individual infection rates Ay entering the expression for R depend on the 
neighbourhood of susceptible nodes as described in Sec. HU 
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Appendix B: Scaling analysis 


In this Appendix, we present the results for the finite-size scaling analysis performed to find the invasion threshold 
in the case of finite- and infinite-range rewiring. 

In order to find the invasion threshold, we exploit the fact that it corresponds to a critical point of the system. 
Following this property, in models with finite-r ang e rewiring, we determined the invasion threshold using finite-size 
scaling for one-dinrensional spanning clusters (291 f32l . 0 I 55 I ]. The relative number of one-dimensional spanning 
clusters, Ni(a,L), exhibits a maximum near the critical value of the inherent transmission rate, a c (see Fig. EJa)), 
and the values of Ni(a, L)L~ e for varying L should collapse onto a single master curve, N\{x), if plotted against 
x = (a — a c )L 1 / V . The exponents v and 6 and critical transmission rate, a c , are found from the scaling collapse (see 
Fig-Elk)). The value of the universal exponent v ~ 1.3±0.1 is consistent with v = 4/3 (62j,[63] known for percolation 
in lattices without rewiring. This is due to the restriction on the maximum rewiring distance in FSS to be much 
smaller and independent of the linear system size, I? max L. 




-6 -4 -2 0 2 4 

(a - a c )L l / u 


FIG. 9: (Color online) (a) The relative number of one-dimensional spanning clusters, Ni(a,L), vs the inherent 
transmission rate, a , for several lattice sizes, L, as indicated in the figure legend, (b) The scaling collapse of 
Ni(a, L)L~ e when plotted against ( a — a c )L X I V with a c = 0.727, v ~ 1.3 ± 0.1 and 9 = —0.0522 for SIR epidemics 
on a square lattice with finite-range rewiring R £ [1,8] with <f> = 0.025 for an exponential form of the rate given by 
Eq. m with synergy parameter, /? = — 1. The different symbols in both panels correspond to different lattice sizes as 
indicated in the legend with each point averaged over 20000 realisations of the epidemics. 


The finite-size scaling used in the case of finite-range rewiring models cannot be applied to infinite-range rewi ring . 
This is due to the existence of an additional length scale related to the distance between the shortcuts [11, 21. 1 641] . 
Therefore, we used two complementary methods for estimating the critical threshold, a c . The first method is based 
on the fact that small-world networks can be described by a mean-field approximation and thus, the relative mass of 
the infinite cluster, A/, depends linearly on (a — a c ) near the critical point. A linear fit [2l| was then used to estimate 
a c (see solid line in Fig. [TUI a)). Alternatively, the critical value of inherent transmission rate can be found from the 
location of the inflection point on the curve for the mass of the infinite cluster, M(a) 0E3 (see vertical dashed 
line through the inflection point in Fig. flOl a) 1. The big error bars seen in Figs. [3][4] are due to the limited world sizes 
available in the small-world simulations. The position of the inflection point gives an upper estimate on the value of 
a c whilst the linear fit provides a lower bound estimate. The restricted system sizes are caused by limited processing 
power available which is required for large worlds with synergistic effects within the continuous-time Kinetic Monte- 
Carlo algorithm. The different tendencies in M(a) with increasing rewiring probability, <f >, for small-world networks 
for different values of the synergy parameter are seen in Fig. [TO] (a) /? < /!*, the system becomes less resilient with 
increasing </>; (b) (3 cs /?*, the mass of the infinite cluster practically does not depend on </>; and (c) /3 > /I*, the system 
becomes more resilient with increasing (f>. 


Appendix C: Technical details in analysis of analytical model 

In this appendix, first we derive an approximate condition for epidemic threshold given by Eq. ([5]) and then give a 
simplified version for it. 

The condition for the epidemic threshold given by Eq. © depends on the bond-percolation threshold and mean 
transmissibility. To leading order in 4> (with (f> < Cl), the bond-percolation threshold for networks with rewiring in a 
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FIG. 10: (Color online) The relative mass of the infinite cluster, M, for an SIR process on a square lattice of size 
N = 201 x 201 with infinite-range rewiring vs the inherent transition rate, a, for different values of synergy 
parameter: (a) /3 = —2, (b) /? = 2 ~ /3* and (c) /3 = 5 and for varying rewiring probability, <f>, as indicated in the 
legends (the same symbols and line styles are used in all panels). An exponential form of the synergistic 
transmission rate given by Eq. m was used. An example of the linear fit used for estimating the critical inherent 
transmission rate is shown by the solid line in (a) for 4> — 0.015 and the vertical dashed line goes through the 

inflection point for the <j> = 0.015 curve. 


finite range is given by, 


T c {(f>) ~ T c 0 (q) - A 9 (R min , -R max )<(> , (Cl) 

where T c o(q) is the bond-percolation threshold for a regular lattice and A q (R m i n , l? max ) > 0 is a mo del-dependent 
constant. 

In order to obtain a linear approximation for the dependence of (T(<^)) on (j>, we first derive the expression for 
degree distribution in a network with rewiring (following Ref. |2C|) and then obtain an approximation for it in the 
case of small <j> <C 1. For example, we consider the case of the square lattice (q = 4) with bonds rewired according 
to the Watts-Strogatz rewiring algorithm as described in Sec. Elne. a SW network. According to this algorithm, 
for each node, Q = q/2 out of q = Q + Q\ bonds attached to this node are rewired with probability </>. This means 
that an arbitrary node has at least Q bonds attached to it. The remaining Q\ = q — Q bonds can be broken by 
rewiring, so that only n\ out of Qi are still attached to the node. The random number n\ is distributed according 
to the binomial distribution, B ni (Q,<j)). In addition to Q + n\ bonds, a random number ri 2 of new bonds can be 
attached to the node as a result of bond rewiring from other nodes. This number is also distributed with binomial 
distribution B n2 (NQ,(f>/N) with NQ being the total number of bonds in the system (the terms ~ 0(./V -1 ) were 
ignored). Therefore, the node degree k = Q + n\ + distribution pk is given by, 

min (k—Q,Q) 

Pk(N , (f>) = B n (Q,(f>)B k _ Q _ n (NQ,(j)/N) , (C2) 

n—0 

if k > Q and p k = 0 otherwise. In the limit of large N 1, the binomial distribution B k -Q- n {NQ, 4>/N) tends to 
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the Poisson one and 


Pk(N,4>) -+p k {4>) 


min 





(■ Q<f>) k - Q ~ n 

( k-Q~n)\ e 


(C3) 


which does not depend on the system size and coincides with the degree distribution obtained for a ring with nodes 
connected to q nearest neighbours [2Cj. The convergence of pk(N , </>) to pk (</>), as demonstrated in Ref. [[20], is rather 
fast and for N ~ 10 3 the numerical data almost perfectly reproduce the limiting distribution given by Eq. (IC3I) . 

The limiting case for small rewiring probabilities, <j> -C 1, follows from Eq. (IC3j) . In this limit, mainly nodes with 
coordination numbers k = q — 1 , q and q + 1 are present in the network and Pk{4>) is given by, 


Pk{4>) — (1 T 2(j)8k,q—l T (-1 • 


(C4) 


The value of (T*,) in the expression for mean transmissibility given by Eq. Q takes into account the synergy effects, 
i.e. that the transmission of infection from a donor to a recipient can occur in the presence of different numbers of 
infected neighbours (excluding the donor) of the recipient (n = 1, 2,..., q — 1) which can affect the transmission rate 
and thus the transmissibility. Within the model of a time-dependent environment with linear transmission rate given 
by Eq. ©, the mean transmissibility with fixed node degree k for /3' > —a/{k — 1) is given by [32J, 


(T fc ) = 1 - e-“ (1 - s k B(p')) 


k -1 


(C5) 


where Sk = Sfc(a, /?', </>) (parameter of the model) is the probability that a neighbour of a recipient node (excluding the 
donor, i.e. one out of k — 1 neighbours) has been infected within the time-interval [—r, r] if the donor became infected 
and started to challenge the recipient at t = 0. The function B(j3') = 1 — (1 — e - ' 3 )//3' increases monotonically from 
0 to 1 with increasing /?' —>• oo. 

Combining Eqs. 0. (ESJ) and (IC4I) gives the desired linear approximation for (T((j))). Introducing this expression 
and Eq. ED into Eq. 0 leads to an approximate condition for the invasion threshold, a = a c (/3', <f>), given by Eq. ([5|) 
which can be transformed to Eq. d. 

The solution of Eq. (Il2bl) significantly depends on the functional form of the infection probabilities s q (a, /?', </> = 0) 
which can be found numerically in the same way as described in Ref. [32j . The dependence of s q (a,(3',(j) = 0) on a 
for q = 3, 4 and 5 and different values of /3' are shown in Fig. |TT| It follows from this figure that s q (a, /3 1 ,<f> = 0) < 1 
for all values of q (at least for the values close to the invasion boundaries). 

As a consequence, the function s q (a, /?' ,0)B(/3') <C 1 and Eq. (I12bl) can be simplified as follows, 


A q -e~ a B(fc) 


T c0 -[l-e-^(l-(q-l)s q (a„pi,0)B(Pi))] =0 
4 {q - l)s 9 (a*,/3(,0) - 2(q - 2)s 9 _i(a*, /?', 0) - 2qs q+1 (a^, /3', 0) + 

ds q (a*,l 3(,0) 


(9-1)- 


d(j> 


= 0 . 


(C6a) 


(C6b) 


Here, the derivative ds q (a*, /?*, 0)/d(/) is typically much smaller than s g (a», /3», 0), as seen from Fig. dd). If the 
expression in the square brackets is positive, Eq. (IC6bl) can be solved, with (a*,/3') lying on the phase boundary for 
<f> = 0 which follows from Eq. (ICGal) . 


Appendix D: Models with added and cut bonds 

In this Appendix, we study analytically and numerically two additional models with just adding and cutting the 
bonds instead of rewiring. 

In a finite-range rewiring model, each bond can be rewired with probability (/>, i.e. the original bond is removed 
from the system and a new bond connecting two nodes within rewiring range is created. Let us modify this model in 
such a way that the original bond is kept in the network in addition to the added bond. This is a model with just 
added bonds [13, S3- F° r small values of the probability of adding a bond, ^add 1, the node degree distribution 
for this model is given by the following equation, 


Pq — (1 4<j6 a dd)^(7,4 T 4(j6 a dd<5g,5 - 


(Dl) 















16 



0.6 
0.5 - 
0.4 - 
£ 0.3 - 
0.2 - 
0.1 - 
0 

0.0 


(b) 


/. «r X x'a 3 ^ 




1 - 

.o L+ **#**’ i 


.*&* 


0.2 


0.4 0.6 

a 


0.8 1.0 




FIG. 11: (Color online) Dependence of (a) S 3 , (b) S 4 and (c) S 5 on inherent rate, a, for different values of the 
synergy parameter, /?', (as marked in the legend in (c)) in square lattices with <f> = 0 in (b) and a small value of 
rewiring probability <fi = 2.5 x 10 -3 and R € [1,4] in (a) and (c). The unlabelled curves show values of /?' varying 
stepwise with unity from /3' = 0 (lowest curve) to f}' = 5 (topmost curve) except for the /3' curve. The figure in (d) 
gives the variation of S 4 (a*,/ 3 () with rewiring probability </>. The open circles in panels (a)-(c) represent the location 
of the critical inherent rate, a c , for each value of /?', whilst the solid circles show the location of a*. 


The corresponding expression for F q is given by 


F q (a,p') = A q -e 


4(1 - s q+1 (a,f3' 1 0)B((3')) q -4(1 - s g (a,/3', 0)U(/3')) 


'U9-1 




t\\q—2 0) 


ok/Add 


— A q ~ e~ a B(p') 


4 {q ~ 1 )s q (a, /?', 0) - 4gs 9 + i(a, /3', 0) + (g - 1) 


ds q (a,/3', 0) 
9(pa,dd 


(D2) 


The dependence of a c ((j)) is given by Eq. © with F q (a,/3') obeying Eq. (ID2I) . The results of its numerical solution 
are shown in Fig. I12f al. It can be seen from this figure that the critical inherent rate decreases with probability </> a dd 
irrespective of the value of /3. Such a monotonic trend is expected and agrees with the fact that Eq. (0), which gives a 
necessary condition for a c to be independent of /3, is not satisfied for any [}. The analytical results are well supported 
by the results of the numerical simulations for 0 a dd 1 and /? < 1 . 

Alternatively, the original rewiring model can be modified in such a way that the original bonds are cut with 
probability <fi cu t but new bonds are not added to the network. This is a model with just removed bonds. The node 
degree distribution for this model is given by the following equation, 


Pq — (1 4(^ cu t)^g i 4 T ^4>cut&q,3 • 


(D3) 
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FIG. 12: (Color online) The dependence of the critical inherent transmission rate on the probabilities (a) </> a dd of 
just adding global bonds and (b) 0 cut of just cutting local bonds on a square lattice. The results were obtained for a 
transmission rate given by Eq. |3]) and adding the bonds within a finite range of R G [1,4]. The symbols refer to the 
simulation results whilst the lines represent the model results given by Eq. m with F q (a, 0') for q = 4 given by 
Eq. (ID2I) in the left panel and Eq. (ID41) in the right panel. The same values of 0' were used in both panels. 


The corresponding expression for F q is given by 


F q {oi, 0') = A g - e“ 


4(1 - *,_!(<*, 0', O)B(0')) q ~ 2 - 4(1 - s q (a, 0', 0 )B{0')) 


'\\q- 1 11 


+ {q-l)B{0'){l~s q {a,0',D)B{0')) q 


,^ q _ 2 ds q {a,0', 0) 


dcj) c 


^A q - e~ a B{0') 


4 {q - 1 )s g (a, 0', 0) - 4 (q - 2)s q -i{a,0', 0 ) + (q - 1) 


ds q (a, 0 1 , 0) 

t^0cut 


■(D4) 


The dependence of a c {<j>) is given by Eq. <0 with F q (a, 0') obeying Eq. (ID4I) . The results of its numerical solution 
are shown in Fig. I12! b). In this case, the critical inherent rate expectedly increases with (f > cut . The numeric results 
are again supportive of the analytics for relatively small values of 0 and show qualitatively the same behaviour for 
0 > 1 . 
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